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Abstract 

We investigate the steady state of a two-dimensional random resistor network 
subjected to two competing biased percolations as a function of the bias strength. 
The properties of the linear and nonlinear regimes are studied by means of 
Monte Carlo simulations. In constant current conditions, a scaling relation is 
found between < R > / < R >q and I/Iq, where < i? > is the average network 
resistance, < R >o the Ohmic resistance and Iq an appropriate threshold value 
for the onset of nonlinearity. A similar scaling relation is found also for the 
relative variance of resistance fluctuations. These results are in good agreement 
with electrical breakdown measurements performed in composite materials. 

The study of electrical and mechanical stability of disordered systems is attracting 
a considerable interest because of its implications on both the material technology § 
and fundamental aspects related to the understanding of the behavior of these systems 
|jl|, II, H, 1^. Indeed, the application of a finite stress (electrical or mechanical) to a 
disordered material generally implies a nonlinear response, which ultimately leads to 
an irreversible breakdown in the high stress limit [jl], |] . Such catastrophic phenomena 
have been successfully studied by using percolation theories [|I], |],||, |],|] . In particular, 
large attention has been devoted to the critical exponents describing the resistance 
and its relative noise in terms of the property of the medium ^ |^, |^. However, 
few attempts have been made so far to describe the behavior of a disordered media 
over the full range of the applied stress ^. On the other hand, such a study 
can provide important information about precursor effects and, more generally, about 
the predictability of breakdown Therefore, a satisfactory understanding of 

breakdown phenomena is demanding Jl], ^. 

Here, our aim is to present a percolative model of sufficient generality to address 
the above issues. Our goal is to provide a theoretical framework to study response and 
fluctuation phenomena under linear and nonlinear regimes in a wide class of disordered 
systems. To this purpose, we analyze the evolution of a random resistor network 
(RRN) in which two competing percolative processes are present, defect generation 
and defect recovery. These processes, which determine the values of the elementary 
network resistances, are driven by the joint effect of an electrical bias and of the heat 
exchange with a thermal bath. The bias is set up by applying a constant current. 
Monte Carlo (MC) simulations are performed to investigate the network evolution in 
the full range of bias values. Indeed, depending on the bias strength, an irreversible 
failure or a steady state of the RRN can be achieved. By focusing on the steady state, 
we analyze the behavior of the average network resistance, < R >, and the properties 
of the resistance fluctuations, as a function of the bias. 



The paper is organized as follows. In Sect. 1 we briefly describe the model used. 
Section 2 presents the results of the MC simulations for the resistance and its fluctu- 
ations. The main conclusions are drawn in Sect. 3 



1 The model 

We study a two-dimensional random resistor network of total resistance R, made of 
Ntot resistors, each of resistance r„, disposed on a square lattice. We take a square 
geometry, NxN, where N determines the linear size of the lattice. For the comparison 
with resistance measurements of thin films, the value of can be related to the ratio 
between the size of the sample and the average size of the grains composing the sample. 
An external bias represented by a constant current /, is applied to the RRN through 
electrical contacts realized by perfectly conducting bars at the left and right hand sides 
of the network. A current in is then flowing through each resistor. The RRN interacts 
with a thermal bath at temperature Tq and the resistances r„ depend linearly on the 
local temperatures, T„, as: 

rn(T„) =ro[l + a(T„-To)] (1) 

In this expression a is the temperature coefficient of the resistance, and T„ is calculated 
by adopting a biased percolation model [TO] as: 



^ neig 

T^ = To + A [r^tl + E {^i^f -rml)] (2) 

Here, Nneig is the number of first neighbours around the nth resistor, the parameter A 
describes the heat coupling of each resistor with the thermal bath and it determines 
the importance of Joule heating effects. The parameter B is taken to be equal to 
3/4 to provide a uniform heating for the case of perfect network. In the initial state 
(corresponding to the perfect network with no heating) all the resistors are identical: 
rn = To- Accordingly, the initial value of the resistance is Rq = ro[N/[N + 1)]. Now, 
we assume that two competing percolative processes act to determine the RRN evo- 
lution. The first process consists of generating fully insulating defects (resistors with 
very high resistance, i.e. broken resistors) with probability Wd^u = gxp[—Ed/ KsTn], 
where is an activation energy characteristic of the defect and Kb the Boltzmann 
constant The second process consists of recovering the insulating defects with 

probability Wji^n = Gxp[—Eji/KBTn], where Eji is an activation energy characteristic 
of this second process. For A 0, Eq. (H) implies that both the processes (defect 
generation and defect recovery) are correlated percolations. Indeed, the probability of 
breaking (recovering) a resistor is higher in the so called "hot spots" of the RRN. On 
the other hand, for A = Eq. (^ yields T„ = Tq, which corresponds to random per- 
colations p, |lT|. The same is true for vanishing small bias values, when Joule heating 



effects are negligible. As a result of the competition between these two percolations 
the RRN reaches a steady state or fails completely with a critical fraction of defect Pc, 
corresponding to the percolation threshold P] . In the first case, the network resistance 
fluctuates around an average value < R >, while in the second case, an irreversible 
breakdown occurs, i.e. R diverges due to the existence of at least one continuous path 
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Figure 1: Average resistance versus tlie external current. Eacli set of data refers to tlie 
values shown in the legend. 



of defects between the upper and lower sides of the network Q . We note that in the 
limit of a vanishing bias (random percolation) and infinite lattices {N — >■ cx)), the 
expression: < Ejy + KbT ln[l + exp^—Ejj/KBT)] provides a sufficient condition 
for the existence of a steady state fill], while it represents only a necessary condition 
in the case of biased percolation. 

The evolution of the RRN is obtained by MC simulations carried out according 
to the following iterative procedure, (i) Starting from the perfect lattice with given 
local currents, the local temperatures are calculated according to Eq. (ii) 
the defects are generated with probability WD,n and the resistances of the unbroken 
resistors are changed as specified by Eq. (|I]); (iii) the currents z„ are calculated by 
solving Kirchhoff's loop equations by the Gauss elimination method and the local 
temperatures are updated; (iv) the defects are recovered with probability Wji^n and 
the temperature dependence of unbroken resistors is again accounted for; (v) R, in 
and Tn are finally calculated and the procedure is iterated from (ii) until one of the 
two following possibilities is achieved. In the first, electrical breakdown occurs (in 
numerical calculations we stop the iteration when R increases over a factor of 10^ with 
respect to its initial value). In the second, the RRN reaches a steady state; in this case 
the iteration runs long enough to allow a fiuctuation analysis to be carried out. Each 
iteration step can be associated with an elementary time step on an appropriate time 
scale (to be calibrated with experiments). In the simulations, as reasonable values of 
the parameters, we have taken: = 75, tq = l(f2), a = 10"^ {K~^), A = 5 x 10^ 
(K/W) and Tq = 300 (K) if not stated otherwise. Furthermore, Ed = 0.17 (eV) 
and values of En ranging from 0.026 (eV) to 0.155 (eV) are used. The values of the 
external bias range from 10~^ < / < 2.8 (A). 
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Figure 2: Log-log plot of the normalized resistances versus the normalized bias. Each set of 
data refers to different values of as reported in the legend of Fig. 1. The insert reports 
the quantity ((i?) — {Rq)) /{Rq) and the dotted line fits the power law y oc with n = 2. 



2 Results 

Figure 1 reports the mean value of the RRN resistance as a function of the external 
current for the values of Er reported in the figure. Each value of < i? > has been 
calculated by consider the time average on a single realization and then averaging over 
20 realizations. The general behavior of each set of data shows a systematic increase 
of the resistance up to the breakdown of the network which occurs for currents greater 
than the breakdown value lb corresponding to the last point reported in the figure. For 
each set of data, at increasing current values, we identify three main transport regimes, 
namely: a linear regime at the lowest currents, a nonlinear regime at intermediate 
current values, and a pre-breakdown regime at the highest currents. The linear regime 
corresponds to a resistance which is independent of current (Ohmic response) and 
extends wider the lower the recovery energy. The nonlinear regime can be described 
by a quadratic dependence of the resistance in an appropriate range of current values 
starting from a threshold value Jq and will be analyzed with details in the following. 
The pre-breakdown regime covers a narrow region just before the breakdown of the 
network which will be better analyzed on the basis of the resistance fluctuations. The 
fact that by increasing the recovery energy the breakdown occurs at lower currents is 
explained by the concomitant growth of the average number of defects, and thus of 
the disorder, inside the network which implies its lower stability against breakdown. 

The similar feature of the data reported in Fig. 1 suggests the possibility that all 
the curves belong to the same universahty class. We define the value of the current /q 
by adopting the criterium of the departure of the resistance from the linear response 
value of 0.05%. Accordingly, Fig. 1 has been redrawn on a normalized current scale 
I / lo as is reported in Fig. 2. This figure, by showing that all curves collapse into a 
single one within an uncertainty of a few percent, clearly confirms the above suggestion. 
Furthermore, to obtain a better insight of the quadratic regime, in the insert we report 
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Figure 3: Log- log plot of the relative variance of resistance fluctuations versus the external 
current. Each set of data refers to the S/j values shown in the legend. 



the relative variation of resistance ((i?) — (i?o)) / (-Ro) as a function of the normalized 
current. Here, we have found that the power law: 

(3) 

is well satisfied over about two decades of current values. The evidence of such a 
quadratic regime is in good agreement with recent experimental results concerning 
the electrical failure of carbon high-density polyethylene composites (HDPE). In par- 
ticular, their failure transition is found to be of the first order which is a feature we 
have checked in the present model. From the study of HDPE samples with different 
concentration of conductive component (carbon black) it was discovered that the ratio 
i?b/i?o between the breakdown resistance i?;, and the Ohmic value Rq is a constant 
independent of the value of the carbon fraction and the geometry of the sample, de- 
pending only on the nature of the conducting component. This result emerges, in our 
context, as a consequence of the independence of the ratio ///q on the initial resis- 
tance [0 and of Eq.(|^). The electrical properties of these composites are thus nicely 
described by a RRN modelling. 

An important source of information concerning the physical properties of the sta- 
tionary regime comes from the study of the resistance fluctuations. To this pur- 
pose we analyze in the following the variance of resistance fluctuations as obtained 
by the simulations. Figure 3 reports the relative variance of resistance fluctuations 
{/\R^)/{Rf = S, with (Ai?2) = (i?2) _ (i?)2^ as a function of the external current for 
the different values of recovery energy reported in the figure. The general behavior of 
different sets are similar to those reported in Fig. 1 for the resistance. Also in this case 
we find a systematic increase of the variance at increasing current with the remarkable 
novelty that now the increase can be of several order of magnitude for the most ordered 
network. By increasing the value of E^, the curves start with higher values, in other 
words, the variance grows with E^. Each curve is initially independent of current. 
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Figure 4: Log-log plot of the relative variance of resistance fluctuations versus pc — p- All 
the values refer to the linear regime. 



which corresponds to the linear regime where the variance is an intrinsic property of 
the RRN. Then, the variance starts increasing with current according with the estab- 
lishment of the non-linear regime driven by the external current and accounted for 
by the bias percolation. Finally, the breakdown of the network is achieved but in a 
much sharp way than for the case of the resistance. We note, that the results at the 
two highest values of the recovery energy (open circles and full squares in the figure) 
are noisy already in the linear regime due to an initial network which is significantly 
disordered and thus less stable than those corresponding to low values of Eji. 

The linear regime, occurring for small currents, when Joule heating effects are 
negligeable, practically corresponds to the competition of two random percolations 

0. We have thus investigated what kind of noise characterizes this random regime. 
Figure 4 reports the relative variance of resistance fluctuations E, as a function of the 
network defectiveness monitored hj pc — p where pc = 0.5 is the critical value of the 
fraction of defects (i.e. the percolation threshold) for a square lattice and p the average 
fraction of defects under stationary conditions. The plotted data have been calculated 
for bias I < Iq and each point refers to the different value of reported in the figure. 
For recovery energies smaller than about 0.103 (eV) we obtain a quasi-perfect lattice, 

1. e. p <C Pc, while for larger values, the fraction of defects become significant and E 
grows following the scaling law: S oc (pc — p)^- Our simulations give k = 3.2 ± 0.1, 
thus confirming previous findings obtained by taking the thermal coupling parameter 



A = |Tl|. We note, that in literature the value 1.12 ± 0.02 has been found for the 
case of a flicker noise spectrum in disordered network subjected to random percolation. 
Therefore, the present noise should be interpreted as coming from the intrinsic source 
of fluctuations and with the property of being much more sensitive to the presence of 
disorder. 

As for the case of the resistance, the similar features of the data reported in Fig. 
3 suggest the possibility that all the curves belong to the same universality class. To 
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Figure 5: Log-log plot of S/Sq — 1 where S/Sq is the relative variance of resistance fluc- 
tuations normalized to the same quantity calculated in the linear regime. The dashed line 
emphasizes the region of validity of the quadratic power law. 



this purpose, Fig. 5 reports on a log-log scale the quantity S/Sq — 1 as a function of 
the normalized current ///q. Here Sq indicates the variance of resistance fluctuations 
for the Ohmic regime. The data are the same of those in Fig. 3 and all the curves are 
found to collapse onto a single one. Looking in detail the data in Fig. 5, we find that 
beside the quadratic regime, higher order terms are present [Q. Indeed, by analogy 



with the case of resistance, the quadratic regime is identified at the lowest current 
(/ < 10 /q) as emphasized by the continuous line. At high currents, (J > 10 Jq) a 
4-th order and eventually a 6-th order regime can be identified. The presence of these 
high order terms can be accounted for by a unique law which expands the variance of 
resistance fluctuations in an even power series of the external current. In particular, 
for large values of the breakdown happens for small /ft, so only the quadratic term 
of the perturbative expansion, [I /I^f , is appreciable. By contrast, for small values of 
Er the breakdown happens for large Ih and terms of order higher than the second are 
detected. We conclude that the noise is more sensitive than resistance to probe the 
pre-breakdown region, since fluctuations are here dramatically amplified with respect 
to negligible deviations of the resistance from the quadratic regime. 



3 Conclusions 

We have studied the steady state of a two-dimensional random resistor network re- 
sulting from the simultaneous evolutions of two competing percolations biased from 
an external current. The two percolations consist in a defect generation process and 
in a defect recovery one. The state of the network is investigated by Monte Carlo sim- 
ulations over the full range of the applied stress. In particular, we have analyzed the 
behavior of the average resistance and of the relative variance of resistance fluctuations 
as a function of the external current and for different values of the recovery activation 



energy. We have found that both these quantities follow a scaling relation in terms of 
///q, i.e. the ratio between the applied current and the current value corresponding 
to the onset of nonlinearity. The scaling exponents are found to be independent of 
the activation energy value, implying that all the failure processes belong to the same 
universality class. In the linear regime, the relative variance of resistance fluctuations 
has been analysed with respect to the average fraction of defects in the network and, 
at increasing disorder, a scaling law as a function of Pc ~ P has been obtained. In the 
nonlinear regime, the relative variance of resistance fluctuations exhibits a quadratic 
dependence on the current while the presence of higher order terms characterizes the 
pre-breakdown region. The analysis of these higher order terms is left as a subject of 
further investigation. 
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